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How does the scaling for the polymer chain in the 
dissipative particle dynamics hold? 
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q ' We performed a series of simulations for a linear polymer chain in a solvent using dissipative particle dynamics 

C/2 | to check the scaling relations for the end-to-end distance, radius of gyration and hydrodynamic radius in three 

dimensions. The polymer chains of up to 80 beads in explicit solvent of various quality are studied. To extract 
the scaling exponent v, the data are analyzed using linear fits, correction-to-scaling forms and analytical fits 
to the histograms of radius of gyration distribution. For certain combinations of the polymer characteristics 
and solvent quality, the correction-to-scaling terms are found to be essential while for the others these are 
negligibly small. In each particular case the final value for the exponent v was chosen according to the best 
least-squares fit. The values of v obtained in this way are found within the interval v = 0.55 4- 0.61 but are 
concentrated mostly around 0.59, which is very close to the best known theoretical result v = 0.588. The 
existence of this interval is attributed both to the peculiarities of the method and to the moderate chain lengths 
being simulated. Within this shortcoming, the polymer chain in this kind of modeling is found to satisfy the 
scaling relations for all three radii being considered 
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q ■ 1. Introduction 

Variety of physical, chemical, and biological phenomena which involve polymers as well as 
' variety of polymer characteristics one is interested in have been naturally reflected in numerous 

theoretical methods and models used for polymer description [ [T] . In particular, it is generally 
recognized by now, that the scaling properties of a long flexible polymer chain immersed in a good 
solvent are perfectly described by the model of self-avoiding walks (SAW) A textbook example 
is given by the mean square end-to-end distance i?uv of a SAW of N steps and that of a linear 
polymer chain, either of which scales in an asymptotic limit TV — > oo as: 



(RIn) ~ N», (1) 

with a universal exponent v. The value of the exponent depends on the space dimension only and 
is the same for all SAW on different three dimensional (3d) lattices. Relation (JIJ is violated for 
the finite N and an approach to the asymptotics is governed by the corrections to scaling. With 
an account of the first correction equation |T]) reads: 

{Rl N )=AN 2 ^l + ^ K + ---y (2) 

Here A, B are non-universal amplitudes and A is an universal correction-to-scaling exponent. 
Theoretical estimates for the exponents follow from the field theoretical renormalization group 
calculations [|3]: 

u(3d) = 0.5882 ± 0.0011, A(3d) = 0.478 ± 0.010. (3) 
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Whereas the SAW on the discrete lattice is a perfect model to study scaling of long flexible 
polymer chains in a good solvent [ |3] , for obvious reasons it fails to describe most of the other 
properties of the chain [ [JJ [5] . A more realistic model would be defined in a continuous space 
and involves both monomer-monomer and monomer-solvent interactions. Brownian dynamics is 
one of the candidates, as it operates on the mesoscale and incorporates the random force acting 
on each monomer which mimics the effect of a solvent [|6|. However, the method lacks correct 
hydrodynamics limit. Similar but more rigorous technique, i.e., the dissipative particles dynamics 
(DPD), was introduced by Hoogerbrugge and Koelman [ [7j and later refined by Espaiiol and 
Warren to satisfy the detailed balance [ [8] . The method has correct hydrodynamics limit [ |9] and 
quickly gained much attention as a powerful technique for mesoscale simulations of various kinds 
of macromolecules [ [6] . 

Therefore, an important issue of analysis is to check how the scaling laws (e.g. of equation ([I]) 
form) hold in the case of the DPD-based simulations. This issue was partly considered in [ ITTJl [TTl 
[12l O [14] . However, some issues remain unclear (see the next section [2] for details). We see a good 
reason to discuss in detail the scaling laws for the polymer chain in a solvent in the DPD method. 
As far as the method was successfully used for the simulations of branched, amphiphilic and other 
complex molecules, the validation of correct scaling laws on such a well-studied system as linear 
chain polymer would be a strong supporting factor for extending the subject to the scaling laws in 
more complex systems (e.g. star-like, branched, dendritic systems, etc.). 

The set-up of our paper is as follows: in the next section [2] we give a brief review of former 
simulation results, section [3] contains the details of our simulations. Our main results are presented 
in section [H section [5] gives conclusions and outlook. 

2. Previous simulations 

A number of various simulation techniques have been used to numerically verify the scaling 
laws of the polymer in solution, namely the molecular dynamics (MD), Monte Carlo (MC), DPD 
and others. Here we shall provide a rather brief account of references relevant to the present study. 

By using the MC method enriched by special techniques, one can move very efficiently through 
the conformational space (see, e.g. [H]). As the result, the method is quite suitable for the simu- 
lations of random walks (RW) and SAW [ |T5l [16] . Li et al. [ [16] used a pivot algorithm and other 
techniques applied to the chain of N = 100 -j- 1000 monomers. The value for the scaling exponent 
obtained is v = 0.588 and it agrees very well with the theoretical estimate ([3]). To study dynamical 
properties, however, one should turn to true dynamical methods. Pierleone and Ryckaert [ [17] 
performed a comprehensive analysis of the scaling laws for both static (radius of gyration, i? g and 
end-to-end distance Rin) and dynamical (diffusion coefficient, D) properties using MD simula- 
tions. One of the important outcomes is that the static properties are unaffected by the finiteness 
of the box size L providing that L/R g > 3. The scaling exponent values v = 0.59 and 0.568 (for the 
Rin and R g , respectively) were obtained. On the contrary, dynamical properties were found to be 
strongly dependent on L. This fact was explained by hydrodynamic interactions of the chain with 
its own images due to the periodic boundary conditions (PBC). If, following Dunweg and Kremer[ 
[T8] one uses the Kirkwood formula, then a correct dynamical scaling is also recovered [117). 

A number of simulations of linear polymers using DPD technique are available. Schlijper et al. 
[ [TO] were the first to study static and dynamic scaling laws for the DPD chain in athermal solvent. 
Later their results for this case were confirmed and extended to include the solvents of variable 
quality by Kong et al. [[UJ. These authors studied R g (among the other properties) and have 
found for the athermal solvent v = 0.52 which is closer to the polymer at the ^-condition. With 
the increase of solvent quality, the value v ss 0.6 close to a SAW exponent was recovered. This 
problem was further addressed by Spenley [[12] who studied the Rin in polymer melt and solution. 
In the latter case, relevant to our study, the exponent v = 0.58 was obtained. The previous result 
v = 0.52 by Kong et al.[[UJ was commented as being influenced by the type of spring potential 
used there. In our opinion, this discrepancy could also be due to insufficient statistics. 

It has also been debated whether the soft repulsion inherent to the DPD method is sufficient 
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to model the self-avoidance of the chains. In particular, Symeonidis et al. [ [13] used additional 
intramolecular interactions (e.g. Lennard- Jones-like repulsion term) to increase the chain self- 
avoidance and experimented with different forms of the bond springs (i.e. FENE, Hookcan and 
WLC) . One of the conclusions made by the authors is, that the presence of FENE potential alone 
ensures a correct value for the exponent v m 0.6 and no extra interactions for the improvement of 
self-avoidance are required. 

Jiang et al. [[14] have carried out an extensive study of the hydrodynamical properties obtained 
via the DPD method. In particular, scaling laws for both R g and Rin for the dilute athermal so- 
lution are relevant to our study. The authors found the values for the scaling exponent v weakly 
dependent on the box size with the averages oiv — 0.58 for the radius of gyration and v = 0.595 for 
the end-to-end distance [1 14] (the box sizes of L = 15-^30 DPD length units were used and a range of 
chain lengths was N = 10 4- 100 DPD beads). The other conclusion given in [[14] is the one already 
mentioned above, claiming that the original soft interactions are sufficient for self-avoidance of the 
chains in the DPD method. Similarly to the MD study of Pierleone and Ryckaert [Hi], the diffusion 
coefficient was found to be highly sensitive to the simulational box size. However, after the appro- 
priate corrections to the system size were introduced, the scaling law for the diffusion coefficient 
D ~ N~ v was obtained. Since, according to Kirkwood theory, the diffusion coefficient D is related 
to the hydrodynamic radius i?h [ HE] , the latter is also shown to obey the correct scaling law. 

Despite the availability of simulational studies that concentrate on static and dynamical scaling 
of the polymer in DPD method, there are still certain points that need to be clarified, namely: 

(i) do all the relevant radii, Rin, R g and i?h scale with the same scaling exponent v! 

(ii) are the corrections to scaling relevant and can these improve the results or make these more 
consistent? 

(iii) in case of a good solvent, will the changes in the polymer-solvent interaction potential lead 
to any detectable drift of the scaling exponent vl 

These three points will be addressed in the current study. 

3. Simulational details 

In our study we closely follow the DPD method as discussed by Groot and Warren [Qj)] . Polymer 
molecule (see, figured]) is modeled as a chain of soft beads connected via harmonic springs. The 
bonding force acting on the i-th bead from its bond neighbour j is: 

P* = -knjhj , (4) 

where = \fij |, Tij = fi — fj, fij = fij/fij and k = 4 is the spring constant. Here and thereafter, 
the length, mass, time and energy (k^T) units are chosen to be equal to unity. The solvent is 




Figure 1. The polymer molecule in a solvent, as modeled in the DPD method. Both monomers 
constituting a polymer chain and solvent molecules are considered as interacting soft beads (grey 
and hollow discs, correspondingly). 
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described in an explicit way, in a form of isolated beads of the same origin as the polymer ones 
(see, figure [1]). Each z-th bead (polymer or solvent one) is subject to a pairwise non-bonded force 
from its j-th counterpart. The force acting on the i-th bead due to the interaction with the j-th 
bead consists of three contributions: 

#V= #g + (5) 
where the conservative F$, dissipative and random PJj; contributions are of the following form: 

pC _ { a(l - n^fij, r y < 1, , s 

ij \ o, m > i, 10 > 

F? = -^(njKhj • %•)%. F% = aw^r^At- 1 / 2 ^. (7) 

Here Vij = Vi — Vj, Hi and Vj being the velocities of the beads, fty is Gaussian random variable: 
{0ij(t)) = 0, (dij(t)dki(t')} = (SikSn + S a S jk )S(t - t'). According to Espanol and Warren [ 8 , the 
dissipative and random force amplitudes are interrelated, w°{rij) = (w R (rjj)) 2 and a 2 — 2j, to 
satisfy the detailed balance. The frequently used analytical form is: 

w°(n j ) = (w*(r ij )f = { £- r «) a ' ;;;•.!; (a) 

We choose the value 7 = 6.75. Parameter a in the conservative force F^j defines maximum repulsion 
between two beads which occurs at complete overlap = 0. Three different a parameters can be 
distinguished, a pp for the polymer-polymer, a ss for the solvent-solvent and a ps for the polymer- 
solvent interactions and all these can be varied. As was shown in [ 119] , these parameters can 
be related to the Flory-Huggins parameter \- To integrate the equations of motion, we used the 
modified velocity- Verlet algorithm by Groot and Warren [ [19] . Let us note that the structural 
properties of the chain in a solvent are defined solely by the conservative forces (O, while the 
dissipative and random forces govern the phase space trajectory. 

In our DPD simulations we considered the following set of chain lengths: N = 5, 7, 10, 14, 20, 
28, 40, 56 and 80 beads in four types of solvent: 

athermal solvent 1, a pp = a ss = 25, a ps = 25, 

athermal solvent 2, a pp = a ss = 33, a ps = 33, , Q s 

good solvent, a pp = a ss = 25, a ps = 20, 

very good solvent, a pp = a ss = 25, a ps = 10. 

The terms "good" and "very good" are indicative rather than exact, since only repulsive forces are 
taken into account in DPD. Hence, "good solvent" means less repulsion for the polymer-solvent 
interaction than that for the polymer-polymer one. The last two cases in equations ([9]) correspond 
to the values £ = —0.2 and £ = —0.6, respectively, in notation of Kong et al. [ fTT] . 

All the simulations were performed in NVT ensemble. At first, the initial chain conformation 
was generated as a RW. Then, it was immersed into a relatively large box filled by solvent and short 
preliminary runs were performed to estimate the approximate R g for the polymer. Thereafter, the 
final run was performed using the cubic box of the linear size L « 5i? g . This size is considered to 
be sufficiently large to eliminate the effects of PBC [[TJ1 [17] ■ The number of solvent beads N s was 
chosen to keep the total reduced density of the solution equal to p* = (N + N s )/V = 3. The time 
step was chosen to be At = 0.04 in reduced units. We performed from 3 • 10 6 (for N = 5) up to 
15 ■ 10 6 (for N — 80) DPD steps for the final run. We did not concentrate on the relaxation times 
in detail, but found the behaviour similar to Jiang et al. [HI], r ~ 0.13 7V 181 who used similar 
DPD parameters for athermal solvent. This brings us to the following estimates for the simulation 
runs: about 5 • 10 4 r for N = 5 and 1.7 • 10 3 r for N = 80 with the estimates for intermediate values 
of N falling in between. The chain conformations were saved after performing every 1000 DPD 
steps and last 3/4 of these were used for the subsequent analysis. As far as all the simulations were 
performed using PBC, the unwrapping of the chain at each time step was performed based on the 
connectivity information. 



4 



Scaling of the DPD polymer chain 



4. Results 

We concentrated our analysis on the three following metric characteristics of the linear polymer. 
The instantaneous squared end-to-end distance at a given time t is defined as: 

R 2 1N (t) = (n-r N ) 2 , (10) 

where fi and rjv are the radius-vectors for the first and last bead of the chain, respectively. The 
squared radius of gyration at time t is defined as: 

1 N 

^W^ivE^-^OM) 2 , ( U ) 

1=1 

where i?cOM is the radius-vector of the polymer center of mass. The inverse hydrodynamic radius 
is defined as follows: 

11^1 

(12) 

r : 



where the sum runs over all different bead pairs All three properties are time averaged (denoted 
hereafter as (■■•)) and then the first two can be square rooted while the third one is inverted to 
provide the estimates for Rin, R g and i?h, respectively. 

One should consider R g and R^ as more useful characteristics as compared to i?uv- First of all, 
Rin is meaningful for the linear polymer or for a certain subchain of the branched molecule (e.g. 
backbone, side chain, star-polymer arm, etc.), while both R g and i?h are universally defined. Apart 
from that, both these properties are directly measurable via small-angle scattering experiments and 
via small-angle dynamic light scattering experiments, respectively (for more discussion on this, see [ 

HO]). 

4.1 . Linear fit and first correction-to-scaling fit 

As was already mentioned above, a number of universal scaling laws have been observed for 
the polymer properties [[5] in the limit of infinite polymer length, e.g. equations |T]). To reproduce 
these properties, one needs to consider rather long chains in the simulations [ I16). Due to the 
polymer coil self-similarity, one would expect that such a coarse-grained approach as DPD should 
be capable of reproducing correct scaling laws for much shorter chain lengths. For the moderate 
polymer lengths, usually the scaling is written in terms of a number of bonds, N — 1, instead of 
the number of monomers N. Also, as is known, the correction-to-scaling terms are important for 
moderate system sizes, as is indicated by the studies of the phase transitions in spin lattice models 
[ EU [22l E3] . Taking into account the first correction-to-scaling term, equation ([2]), one has for the 
(Rj N )- 

{R 2 1N ) = A(N l) 2 " (l + r^—X +■■■)■ (13) 



(N-iy 

We shall work with square rooted values 

Kin = ^rJn) « A^(N - 1)" (l + +•••), (14) 

and keep the linear term only in the correction-to-scaling expansion above. Taking a logarithm of 
both sides of equations (|T4]) , we obtain 



1 { 1 B 

\nTl 1N = -hxA + v]n(N -l)+hx 1 + 



2 v ; V 2 ( N ~ 1) 

A' + v ln(N - 1) + ln (l + (jy f 1)A + • • •) 



(15) 
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where A' and B' are self-explanatory. The same formulae can also be applied for lZ g = «/ (R^) and 

lZ h = l/(l/i?h) (in the latter case one has B' « — B). As a result, the same fit (fPoj) can be used 
for all three quantities under interest, 7^uv, lZ g and IZ^. In the case of linear fit the correction-to- 
scaling amplitude B' is assumed to be zero and one recovers the linear fit form. To perform the 
fits we used the least-squares routine mrqmin from the Numerical recipes book [ 124] . The accuracy 
of the fit was monitored by the cumulative squared deviation per one data point. 



"data 

2 



'"data 2 

— E \K 1N [k] - F(NM) , (16) 

"data ~ L J 

where (N^,TZ 1N [k] ) is the fc-th out of "data data points and F(N) is the fitting function. 

It is generally known that very good statistics for each data point are required to ensure robust 
results for the scaling exponent v. Otherwise the latter depends essentially on the selection of data 
points being used for the fit. This was also found to be the case in the course of the present study. 
Hence, relatively lengthly simulations are performed to improve statistics. To check the consistency 
of the results, we performed a series of fits using a gradually increasing range of intervals in the 
chain length N, e.g. N = 5 -f- 28, N = 5 -j- 40, ■ • • N = 5 -j- 80, which are presented below. 

As was anticipated, both cases of "athermal 1" and "athermal 2" solvents were found to show 
the same averaged properties within the accuracy of the simulations. Therefore, the data for both 
cases have been averaged at each N and marked as "athermal" further in the text. At first glance, 
the data points for both IZin and 1Z g are well fitted by the line in log- log scale (see, figure [5]), but, 
in fact, the best fits were achieved using the first correction-to-scaling form (fl~5)l . We should remark 
here that we were unable to achieve stable results for both 1Zin and lZ g using the four-parameter 
fits, in which all A' , v : B 1 and A in (Tl5|) are fitted simultaneously. Instead, we were forced to fix 
the correction-to-scaling exponent A at certain value and to perform the fit over three remaining 
parameters. We used two choices for A, namely the best theoretical value A = 0.478 [ [3] and 
A = l. 




Figure 2. Data points for IZin (on the left) and lZ g (on the right) and their best fits using the 
form (|15[) obtained from the DPD simulations on a set of chain lengths N = 5 -f- 80 at different 
solvent conditions, see notations on both plots. 



The results presented in table [T] indicate that the correction-to-scaling term for IZin is essential 
in the case of athermal solvent only. In the table, we also show the linear fit performed over the 
largest chain lengths, N = 20 80, for the sake of comparison. One can see that the latter leads to 
almost the same result v — 0.586 as the correction-to-scaling fit v = 0.591, which is indicative of 
the crossover quickly decreasing at N > 20. In the cases of good and very good solvents very small 
amplitudes B' are found and, as the result, linear fits are almost equally good (see, tabled]). For 
1Z g the situation is exactly reversed. The correction to scaling is important for the cases of good 
and very good solvents (see, table [2]). We cannot suggest some definite explanation for this effect, 
except the conclusion that preference should be given to the general form (| 1 5[> . rather than to the 
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linear fit, otherwise the exponent v could be essentially under- or overestimated (see, tables Q] EJ . 



Table 1 . Linear fit and first correction-to-scaling fits to the form ()15[) for IZin with fixed exponent 
A = 0.478 and A = f when using various sets of data points. The best fit for the interval 
N = 5 -I- 80 is underlined. 



fitting 


linear fit 


fit with A = 0.478 


fit with A = 


1 


ranj 


re, N 


V 


x 2 


v \B'\ X 2 


V 


\B'\ 


x 2 


athermal solvent 


5 - 


^ 28 


0.564 


1.7e-5 


0.617 [0.370] 9.4e-7 


0.588 


[0.219] 


1.3e-6 


5 - 


^ 40 


0.569 


3.1e-5 


0.621 [0.399] 9.0e-7 


0.592 


[0.252] 


1.8e-6 


5 - 


^ 56 


0.573 


5.5e-5 


0.627 [0.447] 1.3e-6 


0.596 


[0.295] 


3.6e-6 


5 - 


~ 80 


0.574 


5.2e-5 


0.611 [0.318] 8.2e-6 


0.591 


[0.248] 


6.9e-6 


20 


80 


0.586 


8.6e-6 










good solvent 


5 - 


^ 28 


0.607 


8.6e-6 


0.591 [-0.099] 6.4e-6 


0.600 


[-0.066] 


6.9e-6 


5 - 


^ 40 


0.607 


7.7e-6 


0.595 [-0.076] 5.7e-6 


0.601 


[-0.056] 


6.0e-6 


5 - 


^ 56 


0.607 


7.3e-6 


0.604 [-0.026] 7.0e-6 


0.605 


[-0.023] 


7.0e-6 


5 - 


~ 80 


0.605 


1.7e-5 


0.594 [-0.087] l.le-5 


0.600 


[-0.073] 


1.2e-5 


20 


~ 80 


0.601 


1.5e-5 










very good solvent 


5 


-r 28 


0.614 


6.3e-7 


0.610 [-0.021] 5.4e-7 


0.612 


[-0.019] 


5.0e-7 


5 


-r- 40 


0.613 


1.7e-6 


0.605 [-0.052] 8.3e-7 


0.609 


[-0.041] 


8.4e-7 


5 


-j- 56 


0.613 


2.1e-6 


0.612 [-0.009] 2.1e-6 


0.612 


[-0.011] 


2.1e-6 


5 


-j- 80 


0.613 


2.0e-6 


0.611 [-0.012] 1.9e-6 


0.612 


[-0.012] 


1.8e-6 


20 


-T- 80 


0.613 


2.8e-6 











The four-parameter fit, mentioned above, does not work for either TvLuv or TZ g but is supposed 
to perform better for IZh, where essentially larger correction to scaling is to be expected [ [20], 
Indeed, this is found to be the case in our simulations. First of all, this can be noticed from the 
arrangement of the data points, see figure [3] The results based on numerical fits are presented in 
table [3] where both the cases of fixed A = 0.478 and of the four-parameter fits are shown. The 
latter give essentially better results in terms of fitting accuracy and provide independent estimates 
for both exponents, v w 0.56 -j- 0.58 and A w 0.65 -j- 0.70. As was demonstrated by Dunweg et al. 
[ |2"U] there are two comparable correction-to-scaling terms for the 7?.h, one is proportional to N~ A 
and another, the so-called "analytic" correction, is proportional to TV"' 1 " 1 '). Therefore, if one uses 
the form l|15p with a single correction term, then a rather effective exponent A will be obtained. 
This is how we should interpret the value A obtained in our study. 

4.2. The fits based on the distribution of "R e values 



The averaging of the metric properties of interest, e.g. 7Z g = w (i? 2 ), throughout this study was 

performed as the arithmetic mean over discrete set of time points. One can also be interested in 
the distribution p(R g ) of the R g values themselves and, for instance, in examination of the scaling 
law for the most probable values with the chain length N. The following analytical form has been 
postulated by Lhuillier in 3d [ |2"5] 

p(R g ) ~ exp [—Ai(N u / R g ) 3a - A 2 {RJN V ) 5 ] , (17) 

where a — 1/(3^ — 1), 5 = 1/(1 — v) and A\, A 2 are constants. Subsequently, this has been verified 
in a number of lattice MC simulations [[551 HI] • O ne can rewrite this distribution in terms of i? 2 , 
using the expressions for a and S and incorporating the iV-dependent factors into the constants, 
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Table 2. Linear fit and first correction-to-scaling fits to the form ()15|) for 1Z S with fixed exponent 
A = 0.478 and A = 1 when using various sets of data points. The best fit for the interval 
N = 5 -7- 80 is underlined. 



fitting 


linear fit 




fit with A = 0.478 




fit with A = 


1 




ran 


ge, N 


V 


x 2 




v [B>] X 2 




V 


[B'\ 


x 2 




athermal solvent 


5 


~ 28 


0.601 


2.7e- 


-6 


0.584 [-0.106] 2.0e- 


-7 


0.592 


[-0.080] 


2.5e- 


-7 


5 


+ 40 


0.600 


3.1e- 


-6 


0.587 [-0.090] 3.0e- 


-7 


0.593 


[-0.073] 


2.6e- 


-7 


5 


-r 56 


0.600 


2.8e- 


-6 


0.593 [-0.052] 1.5e- 


-6 


0.597 


[-0.047] 


1.3e- 


-6 


5 


~ 80 


0.599 


8.0e- 


-6 


0.587 [-0.087] 2.9e- 


-6 


0.593 


[-0.077] 


3.1e- 


-6 


20 


80 


0.596 


6.0e- 


-6 














good solvent 


5 


+ 28 


0.641 


4.9e- 


-5 


0.573 [-0.368] 2.7e- 


-6 


0.605 


[-0.318] 


4.2e- 


-6 


5 


+ 40 


0.638 


5.3e- 


-5 


0.586 [-0.314] 5.4e- 


-6 


0.610 


[-0.283] 


5.0e- 


-6 


5 


+ 56 


0.636 


5.1e- 


-5 


0.597 [-0.260] 9.4e- 


-6 


0.615 


[-0.244] 


7.3e- 


-6 


5 


~ 80 


0.631 


9.7e- 


-5 


0.588 [-0.306] 1.4e- 


-5 


0.609 


[-0.299] 


1.4e- 


-5 


20 


80 


0.618 


2.8e- 


-5 














very good solvent 


5 


+ 28 


0.651 


6.9e- 


-6 


0.626 [-0.144] 1.9e- 


-6 


0.639 


[-0.109] 


2.2e- 


-6 


5 


~ 40 


0.647 


1.6e- 


-5 


0.618 [-0.185] 2.6e- 


-6 


0.633 


[-0.151] 


3.6e- 


-6 


5 


~ 56 


0.645 


2.4e- 


-5 


0.616 [-0.194] 2.4e- 


-6 


0.631 


[-0.168] 


3.7e- 


-6 


5 


+ 80 


0.642 


3.4e- 


-5 


0.615 [-0.201] 2.2e- 


-6 


0.629 


[-0.184] 


3.9e- 


-6 


20 


80 


0.632 


2.3e- 


-6 















Table 3. Linear fit and first correction-to-scaling fits to the form (|15[) for TZh with fixed exponent 
A = 0.478 and from the four-parameter fits when using various sets of data points. The best fit 
for the interval N = 5 4- 80 is underlined. 



fitting 
range, N 


fit with A = 0.478 
v \B'\ X 2 


four-parameter fit 
v A [B'\ x 2 


athermal solvent 


5 - 
5 - 
5 - 
5 - 


- 28 

- 40 

- 56 

- 80 


0.596 [3.008] 4.9e-5 
0.611 [3.275] 6.9e-5 
0.625 [3.555] 9.5e-5 
0.635 [3.786] 1.0e-4 


0.497 0.836 [1.324] 1.9e-7 
0.534 0.732 [1.683] 4.7e-7 
0.572 0.666 [2.195] 1.0e-6 
0.552 0.699 [1.912] 1.2e-6 


good solvent 


5 - 
5 - 
5 - 
5 - 


- 28 

- 40 

- 56 

- 80 


0.623 [2.953] 2.0e-5 
0.639 [3.274] 3.1e-5 
0.654 [3.592] 5.1e-5 
0.658 [3.710] 4.7e-5 


0.476 1.017 [0.925] 4.8e-7 
0.563 0.654 [1.636] 2.3e-6 
0.641 0.577 [2.993] 3.0e-6 
0.565 0.664 [1.663] 5.5e-6 


very good solvent 


5 - 
5 - 
5 - 
5 - 


- 28 

- 40 

- 56 

- 80 


0.625 [3.052] 5.9e-5 
0.639 [3.284] 7.8e-5 
0.651 [3.533] 9.7e-5 
0.662 [3.779] l.le-4 


0.652 0.643 [3.159] 1.5e-7 
0.580 0.703 [1.925] 6.1e-7 
0.583 0.698 [1.963] 5.4e-7 
0.577 0.710 [1.885] 5.1e-7 
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1.5 



In R h 




• very good 

1.0 - * good 

o athermal 

0.5 - 



0.0 



1.0 1.5 2.0 2.5 3.0 3.5 4.0 ln(JV-l) 

Figure 3. Data points for IZh and their best fits using the form (|15|l obtained from the DPD 
simulations on a set of chain lengths N — 5 4- 80 at different solvent conditions, see notations 
on the plot. 



and the result reads: 



p(i?2) = Cexp 



This distribution has a maximum at 



2imax 



a i 



a 2 [Rl}^^ 



3ai(l — v) 



a 2 (3f - 1) 



(18) 



(19) 



which can be square rooted to obtain an estimate for the most probable value of the radius of 
gyration ft g max = y , [^] mas . We prefer to work with the distribution of R* $J5§ and then to 

square root the maximum to match the averaging procedure for simple averaging, lZ g — <J {Rg) ■ 



From the point of view of numerical fitting, the form (JT8J) has four parameters: one exponent 
v and three coefficients, a\, et 2 and C. Therefore, one can have a first estimate for the exponent 
v right from the fit of the data to the form (JTHJ) ■ The second estimate for v can be obtained from 
the anticipated scaling of 7^ g max according to the form (fT5")) . 




Figure 4. Mapping of the distribution of radius of gyration values p(Rg) obtained from the 
DPD simulations in athermal solvent, onto Lhuillier form (|18p . a few indicated chain lengths are 
shown only for the sake of clarity. 

The distributions p(Rg) were built in a form of histograms and were found to map very well on 
the Lhuillier form (fT5)) . see samples in figureS] obtained for the case of athermal solvent. These fits, 
as was mentioned above, provide both the estimate for lZ g ma * and an independent estimate for 
the exponent v. The latter was found to be scattered in the interval v = 0.53 4- 0.58 depending on 
N and solvent type. We do not expect high accuracy from such a fitting, as far as dependence on 
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v in equations (fl~8|) is rather complex. Nevertheless, the values of v found are very reasonable and 
indicate a self-consistency of this numerical approach. The next step is to use the maxima positions 
7?.g max of the distributions (fl~8|) . found via equations (fT9)l , and to find a fit to the scaling law (fT5|) 
similarly to the analysis of lZ g data in the previous subsection. The data points alongside with the 
best fits are shown in figure [5] and the results for the exponent v, correction-to-scaling amplitude 
B' and fitting error x 2 are presented in table [4] We would like to stress that, contrary to the results 
given in table [2J no dependence on solvent quality is observed in table 2] and the same value v = 
0.582 — 0.585 is found essentially in all three cases. This might indicate that the effect of small drift 
of v observed in a very good solvent in table [2] could be an artefact of the numerical method used. 

2.0 
1.5 
1.0 
0.5 
0.0 
-0.5 

1.0 1.5 2.0 2.5 3.0 3.5 4.0 ln(JV-l) 

Figure 5. Data points for 7?,g max and their best fits using the form (TT5J) obtained from the DPD 
simulations on a set of chain lengths N = 5 ~ 80 at different solvent conditions, see notations 
on the plot. 



Table 4. Linear fit and first correction-to-scaling fits to the form ()15|l for i?™ ax with fixed 
exponent A = 0.478 and A = 1 when using various sets of data points. The best fit for the 
interval N = 5 -f- 80 is underlined. 



fitting 


linear fit 




fit with A = 0.478 




fit with A = 


1 




range, TV 


V 


x 2 




V 


[B>] x 2 




V 


\B'\ 


x 2 














at her 


mal solvent 












5 - 


- 28 


0.602 


6.4e- 


-5 


0.539 


[-0.345] 2.6e- 


-5 


0.568 


[-0.306] 


2.4e- 


-5 


5 - 


- 40 


0.601 


5.6e- 


-5 


0.562 


[-0.242] 3.2e- 


-5 


0.579 


[-0.222] 


2.8e- 


-5 


5 - 


- 56 


0.601 


4.9e- 


-5 


0.576 


[-0.169] 3.4e- 


-5 


0.587 


[-0.165] 


3.0e- 


-5 


5 - 


- 80 


0.598 


5.6e- 


-5 


0.574 


[-0.183] 3.1e- 


-5 


0.585 


[-0.181] 


2.7e- 


-5 


20 - 


■- 80 


0.598 


2.3e- 


-5 
















good solvent 


5 - 


- 28 


0.646 


4.6e- 


-5 


0.594 


[-0.291] 2.1e- 


-5 


0.618 


[-0.247] 


2.0e- 


-5 


5 - 


- 40 


0.641 


6.6e- 


-5 


0.588 


[-0.316] 1.8e- 


-5 


0.613 


[-0.283] 


1.9e- 


-5 


5 - 


- 56 


0.640 


5.9e- 


-5 


0.606 


[-0.231] 2.7e- 


-5 


0.622 


[-0.217] 


2.5e- 


-5 


5 - 


- 80 


0.633 


1.8e- 


-4 


0.582 


[-0.352] 6.2e- 


-5 


0.607 


[-0.344] 


6.9e- 


-5 


20 - 


■- 80 


0.609 


9.4e- 


-5 


























very { 


;ood solvent 












5 - 


- 28 


0.646 


2.6e- 


-4 


0.517 


[-0.626] 6.0e- 


-5 


0.571 


[-0.634] 


5.5e- 


-5 


5 - 


- 40 


0.636 


3.4e- 


-4 


0.522 


[-0.610] 5.2e- 


-5 


0.570 


[-0.637] 


4.7e- 


-5 


5 - 


- 56 


0.628 


3.8e- 


-4 


0.529 


[-0.582] 4.9e- 


-5 


0.572 


[-0.624] 


4.2e- 


-5 


5 - 


- 80 


0.625 


3.7e- 


-4 


0.549 


[-0.503] 7.8e- 


-5 


0.582 


[-0.551] 


5.9e- 


-5 


20 - 


■- 80 


0.607 


3.9e- 


-5 

















The cumulative outcome of our analysis for the scaling exponent v which governs the scaling 
law ([2]) is presented in figure [6l Here we show a histogram of the values of v obtained for 7?.iAr, 




10 



Scaling of the DPD polymer chain 



TZ g and IZ^ of a polymer chain in different solvents considered in the current study, namely, in 
athermal, good and very good solvents. To make this histogram, the best fits (underlined results in 
tables [T] [4j) have been used. One can see from the histogram that the value of v is scattered in an 
interval v — 0.55-^0.61 but in fact is found predominantly in much narrower interval v = 0.58-^0.60 
centered around the best theoretical estimate v = 0.5882 ± 0.0011 [ 3 . This allows us to conclude 
that within the accuracy of our analysis, the scaling law ([2]) holds for all polymer characteristics 
and in all solvents considered in this study. Our conservative estimate for the scaling exponent 
therefore is 

v » 0.59 ± 0.01. (20) 




Figure 6. A histogram of values of v obtained for different characteristic radii of a polymer 
chain in different solvents (|9}. The data are collected according to best fits results (underlined 
in tables [TJ-[4j| . The value of v is found predominantly in the interval v = 0.55 ~ 0.61 and 
concentrate around the best theoretical estimate v — 0.5882 ± 0.0011 [[3]. 



5. Conclusions and outlook 

In this study we performed DPD simulations of the polymer chain in a solvent of various quality 
with the aim to reexamine how well the scaling laws hold for various polymer characteristics. Chains 
of up to 80 beads were considered and simulation runs from 1600 up to 50000 relaxation times 
were performed. Three metric properties were considered, end-to-end distance, radius of gyration 
and hydrodynamic radius. For the analysis we used the linear fits, first correction-to-scaling fits 
and fitting the distribution of the radius of gyration. As an outcome, within the accuracy of the 
simulations and of the data analysis technique, the following conclusions can be drawn: 

(i) All three metric properties obey the scaling law with very close values for v found within 
the interval 0.55 0.61 (depending on the method of analysis). Most values are concentrated 
around the average v — 0.59 being very close to the theoretical estimate v — 0.5882 ± 0.0011 
[[3], see, figure El 

(ii) Corrections to scaling are found to be important in all combinations except a few prop- 
erty/solvent ones and thus the conclusion (i) is valid only when the correction-to-scaling 
terms are taken into account. 

(iii) No or very small (up to 4%) drift of the scaling exponent v is observed when changing the 
solvent quality (remaining, however, in a good solvent regime), but this effect depends on a 
numerical technique being used. 

The method of analysis can be extended to the study of the scaling laws in more complex 
molecular architectures, e.g. star-like polymers, branched and hyperbranched molecules (including 
amphiphiles). 
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